import random, numpy, pickle, math

try:
    factorial = math.factorial
except AttributeError:
    try:
        import scipy
        factorial = scipy.factorial
    except ImportError:
        factorial = lambda n: reduce(lambda a,b: a*(b+1), range(n), 1)


# Testing 

#spf = open('stabpickle','r')
#stabs = pickle.load(spf)
#spf.close 

def fixtab():
    badIN = '#lllllll#'
    alts = altERCs(badIN)
    for I,V in stabs.items():
        # print 'Input:',I
        if I == badIN:
            print 'FIXING',I
            for O,E in V['ercs'].items():
                # print 'O:',O
                Ox = O[1:-1]
                if Ox in alts: pass
                else: print Ox,alts
                V['ercs'][O] = alts[Ox]
    


def testab():
    bad = []
    for I,V in stabs.items():
        print '------------------------------------\nInput:',I
        for O,E in V['ercs'].items():
            good = eBitsConsistent(E,10)
            print '       ',O,good
            if not good: 
                bad.append((I,O,E))
    return bad
    

def search():
    info = {}
    for I,V in stabs.items():
        sc = V['surfaceCandidates']
        dc = V['ercs']
        LnS = len(sc)
        LnD = len(dc)
        if LnS != 9 or LnD > 12: continue  
        # print '------------------------------------\nInput:',I, LnS,LnD
        for (o,cx) in sc.items():
            N = len(cx)
            E = removeSimpleEnts(V['ercs'][cx[0]],10)
            if N > 1: 
                for surf in cx[1:]:
                    E = eBitDisjPair(E,V['ercs'][surf],10)
                # print E
                info[(I,o)] = (len(cx),E)
    # ofInterest = ['#lhh#', '#llh#', '#lhl#', '#lll#']
    ofInterest = ['#lhh#', '#lll#']
    lasDs = 0
    for ((DS,SS),(N,E)) in sorted(info.items()):
        if DS not in ofInterest: continue
        if DS != lasDs:
            print '----------------------------------------------'
            lasDs = DS
        print DS,SS,'vs.',
        xsp = 0
        for DS2 in ofInterest:
            if DS == DS2: continue
            T2 = stabs[DS2]
            E2x = T2['ercs'].values()
            compat = 0
            for E2 in E2x:
                if eBitsConsistent(E|E2,10): compat += 1
            print ' '*xsp,DS2,'yields',compat,'of',len(E2x)
            xsp = 19
    print '================================================='
    for ((DS,SS),(N,E)) in sorted(info.items()):
        if DS not in ofInterest: continue
        if DS != lasDs:
            print '----------------------------------------------'
            lasDs = DS
        for DS2 in ofInterest:
            if DS == DS2: continue
            T2 = stabs[DS2]
            E2x = T2['ercs'].items()
            print DS,SS,'in the tableau for',DS2,'yields:'
            for (C,E2) in E2x:
                if eBitsConsistent(E|E2,10): print ' '*20,'| | live',C
                else: print ' '*20,'|X| DEAD',C
    


    
    
#------------------------------------------------------------------------
# Bit-ERCs are integers whose even bits denote Ws and odd bits denote Ls. 
# E.g., W:10, L:01, E:00, WLEW:10010010=290.  Bitwise or '|' does fusion.
#       W:2   L:1   E:0
#
# But,  W:00, L:11, E:10, is obtained by flipping all the even bits
#       W:0   L:3   E:2   and now entailment is equivalent to '>'.
#
# A -> B if A^m2 > B^m2   Masks: m1 = 4**k/3 = '...010101' 
#                                m2 = m01<<1 = '...101010'
#------------------------------------------------------------------------

def ercBit(P,Q):
    R = 0
    for (ix,p) in enumerate(P):
        R <<= 2
        if p>Q[ix]: R |= 1
        elif p<Q[ix]: R |= 2
    return R


m_Dict = {}
def getM1(k): 
    try: return m_Dict[('m1',k)]
    except KeyError: m_Dict[('m1',k)] = int(4**k/3)
    return m_Dict[('m1',k)]
    
def getM2(k):
    try: return m_Dict[('m2',k)]
    except KeyError: m_Dict[('m2',k)] = getM1(k) << 1
    return m_Dict[('m2',k)]
    

#------------------------------------------------------------------------
# Translation to and from other representations
#------------------------------------------------------------------------

def ercStr2bit(ercStr):     # strings to bits
    R = 0
    for p in ercStr:
        R <<= 2
        if p=='l': R |= 1
        elif p=='w': R |= 2
    return R

def ercBit2str(eBit,k=None):   # if k is not specified initial E's are dropped 
    E = ''
    while eBit:
        if eBit&1: E = 'l' + E 
        elif eBit&2: E = 'w' + E  
        else: E = 'e' + E 
        eBit >>= 2
    if k: E = ((k - len(E))*'e') + E
    return E

#------------------------------------------------------------------------
# Random lists of violation vectors and eBits for testing.
#------------------------------------------------------------------------

def ranDat(k,n):
    R = []
    for row in range(n):
        r = []
        for con in range(k): r.append(random.choice([0,1,2]))
        R.append(tuple(r))
    return R

def ranBits(k,n):
    D = ranDat(k,n)
    return [ercBit(D[0],Q) for Q in D]


def ranBitsConsistent(k,n):
    D = ranDat(k,n)
    EB = [ercBit(D[0],Q) for Q in D]
    if eBitsConsistent(EB,k): return EB
    else: return ranBitsConsistent(k,n) 

#------------------------------------------------------------------------
# Consistency check for a set of eBits
#------------------------------------------------------------------------

# def eBitsConsistent(eBits,k):
#     print 'eBits:', eBits
#     M = numpy.ma.array(eBits,int)
#     #if m1 == None: 
#     #    m1 = 4**k/3  # phil is a length-k '...010101' mask
#     # m1 = getM1(k)
#     # m2 = getM2(k)
#     m1 = 4**k/3
#     m2 = m1<<1
#     for x in xrange(100):
#         print 'M before:', M
#         for y in M:
#             if y:
#                 print numpy.binary_repr(int(y))
#         Fuse = numpy.ma.bitwise_or.reduce(M)
#         print 'Fuse:', Fuse,numpy.binary_repr(int(Fuse))
#         Lose = Fuse&m1
#         print 'Lose:', Lose, numpy.binary_repr(int(Lose))
#         if not Lose: return True
#         Winz = Fuse&(m2)
#         print 'Winz:', Winz, numpy.binary_repr(int(Winz))
#         Next = (Winz^(Lose<<1))&Winz
#         print 'Next:', Next, numpy.binary_repr(int(Next))
#         if not Next: return False  
#         M = numpy.ma.masked_where(M&Next, M)
#         print 'M after:', M

def opad(b,k):
    bnr = numpy.binary_repr(b)
    if len(bnr) < 2*k: return 'O'*(2*k-len(bnr))+bnr
    else: return bnr
    
def eBitsConsistent(eBits,k):
    M = numpy.ma.masked_array(list(eBits))
    verbose = False
    m1 = getM1(k)
    m2 = getM2(k)
    while M.any():
        Fuse = numpy.ma.bitwise_or.reduce(M)
        if Fuse == '--': return True
        if verbose: 
            print '-----------------\n', M, '\n-----------------'
            print ' Fused = binary',opad(Fuse,k), 
            print 'erc',ercBit2str(Fuse,k),
            print 'int:',Fuse
        if not Fuse:
            return True
        Lose = Fuse & m1
        if not Lose: 
            return True
        WW = Fuse & m2
        Winz = (WW^(Lose<<1))&WW
        if verbose: 
            print 'L-vals = binary',opad(Lose,k), 
            print 'erc',ercBit2str(Lose,k),
            print 'int:',Lose
            print 'W-vals = binary',opad(Winz,k), 
            print 'erc',ercBit2str(Winz,k),
            print 'int:',Winz
        if not Winz: return False  
        M = numpy.ma.masked_where(M&Winz, M)

    return True


def tester(k,n):
    D = ranDat(k,n)
    for row in D: print row
    print 'Ercs for first row vs. each row (including itself):'
    EE = [ercBit(D[0],Q) for Q in D]
    ES = [ercBit2str(eBit,k) for eBit in EE]
    print EE
    print ES
    for e in EE: print ercBit2str(e,k),
    print 
    
    Tv1 = eBitsConsistent(EE,k)
    Tv2 = consistent(ES)
    
    if Tv1 != Tv2: 
        print '\noh Shit!'
        print Tv1,EE,Tv2,ES
        return True
    else: 
        return False
    
    #  print '\nRow 1 is a contender.'
    # else: print '\nRow 1 is harmonically bounded.'
    # if consistent(ES): print 'Really it is!'
    # else: 'Really it is not!' 


def testy(n):
    for x in range(n):
        if tester(6,6): break
    print 'Finished all',n

#-----------------------------------------------------------------------
# Consistency check for a bunch of violation vectors at the same time. 
#-----------------------------------------------------------------------

def matCon(vioVex,k):
    '''Takes a list of violation profiles and returns the contenders.'''
    phil = 4**k/3 # phil is a length-k '...010101' mask
    M = numpy.ma.array([[ercBit(P,Q) for P in vioVex] for Q in vioVex])
    while 1: 
        fused = numpy.ma.bitwise_or.reduce(M,axis=0)
        L = fused&phil
        W = fused&phil<<1
        S = (W^L<<1)&W
        print '-------------\n%s\nF: %s\nL: %s\nS: %s' % (M, fused, L, S)
        if any(S): M = numpy.ma.masked_where(M&S, M)
        else: return [v for (i,v) in enumerate(vioVex) if L[i]==0] 



def testerM(k,n):
    D = ranDat(k,n)
    for row in D: print row
    return matCon(D,k)

#-----------------------------------------------------------------------
# Let's get some ERC-hulls
#-----------------------------------------------------------------------

def eBitOr(eBit1,eBit2,k):     # ERC disjunction equals union of
    m2 = getM2(k)
    L = eBit1&eBit2             # the Ws & intersection of the Ls
    W = (eBit1&m2)|(eBit2&m2)   # -- m2 is the ...101010 bit mask    
    return L|W


def eBitDisj(ebSets, k):
    E = ebSets[0]
    for E2 in ebSets[1:]: 
        E = eBitDisjPair(E,E2,k)
    return E 

def eBitDisjPair(ebSet1, ebSet2, k): 
    m1 = getM1(k)
    S = set([])
    for e1 in ebSet1:               
        for e2 in ebSet2:         
            E = eBitOr(e1,e2,k)
            if E&m1: S.add(E)     # skip disjunctions lacking L
    R = removeSimpleEnts(S,k)
    return R


#-----------------------------------------------------------------------
# Entailment among Bit ERCs
#-----------------------------------------------------------------------

def removeSimpleEnts(SS,k):
    S, R = set(SS), set([])
    while S: 
        e = S.pop()
        if unentailed(e,S,k) and unentailed(e,R,k): R.add(e)
    return frozenset(R)

def unentailed(e,Eset,k):
    '''Returns True if no erc in Eset entails e.'''
    m2 = getM2(k)
    for ee in Eset:
        if eBitEntails(ee,e,k): return False
    return True

def eBitEntails(A,B,k):
    '''Returns True if eBit A entails eBit B.'''
  # if (A has unique W's) or (B has unique L's) then not A->B
    if ((A^B)&A&getM2(k)) or ((A^B)&B&getM1(k)): return False
    else: return True
        

def tesDisj(k,n):
    A,B,C = ranBitsConsistent(k,n), ranBitsConsistent(k,n), ranBitsConsistent(k,n)
    for a in A: print a, ercBit2str(a,k),
    print 
    for b in B: print b, ercBit2str(b,k),
    print 
    for c in C: print c, ercBit2str(c,k),
    print 
    R = eBitDisj([A,B,C],k)
    for r in R: print r, ercBit2str(r,k),
    print 
    return R


# Volume computation, what's up.

def countWL(eBit):
    cW, cL = 0, 0
    while eBit:
        if eBit & 1:
            cL += 1
        elif eBit & 2:
            cW += 1
        eBit >>= 2
    return cW, cL

def shift(eBit, col, k):
    ebShift = eBit
    nShifts = 2*(k - col - 1)
    return ebShift >> nShifts

def withoutCol(eBit, col, k):
    ri = k - col
    nShifts = 2*ri
    left = ((eBit >> nShifts) << nShifts) >> 2
    rightmask = 3
    while ri:
        rightmask = (rightmask << 2) | 3
        ri -= 1
    right = eBit & rightmask
    
    print 'ebit: ', bin(eBit)
    print 'left: ', bin(left)
    print 'mask: ', bin(rightmask)
    print 'right:', bin(right)
    return left | right


def eBitsRVol(ercs, k):
    if not eBitsConsistent(ercs, k):
        return 0.

    def __rVol(eBits, ercLen):
        volume = 0.
        if len(eBits) > 1:
            for i in xrange(ercLen):
                colHasL = False
                colHasE = False
                for eb in eBits:
                    ebShift = shift(eb, i, ercLen)
                    if ebShift & 1:
                        colHasL = True
                        break
                    if ebShift & 2:
                        colHasE = True

                if colHasL: pass
                elif colHasE:
                    recurse = __rVol([withoutCol(eb, i, ercLen) for eb in eBits
                        if not (shift(eb, i, ercLen) & 2)], 
                            ercLen-1)
                    volume += recurse
                else:
                    volume += factorial(ercLen-1)
        else:
            cW, cL = countWL(eBits[0])
            ploar = float(cW+cL)
            kfac = factorial(ercLen)
            if polar: volume = (cW/polar)*kfac
            else: volume = kfac
        return volume

    return __rVol(ercs, k)


#----------------------------------------------------------------------
# Testing against the old representation
#----------------------------------------------------------------------

def randomData(i,j):
    '''Random i-by-j matrix of integers 0-9.'''
    X,Y = range(i), range(j)
    R = [[random.randint(0,9) for x in X] for y in Y]
    return R
            
def ERCs(n,M):
    """
    An ERC-set is a set of statements in a 3-valued logic expressed with
    n-tupples from {w,l,e}^n -which are represented as strings here. The
    meaning of an ERC is that at least one of the constraints with a 'w' 
    outranks all of the constraitns with 'l'. E.g. 'wlwel' means that C1
    or C3 outranks both C2 and C5. This function retunrs the ERC-set for   
    the n-th row of integer matrix M (where n for the first row is 0).
    """
    R = []
    for rowIndex in range(len(M)):
        if rowIndex != n: 
            e = ''
            for (v1,v2) in zip(M[n],M[rowIndex]): 
                if v1>v2:  e+='l'
                if v1==v2: e+='e'
                if v1<v2:  e+='w'
            R += [e]
    return R   

def fuse(e1,e2):
    R = ''
    for (a,b) in zip(e1,e2):
        if a =='l' or b =='l': R+='l'
        elif a=='e' and b=='e': R+='e'
        else: R+= 'w'
    return R

def fusion(eset):
    '''Fusion of an entire ERC set.'''
    return reduce(fuse,eset)

def consistent(eSet):
    '''RDC-style consistency check.'''
    eset = list(set(eSet))
    if len(eset) > 1: 
        fused = fusion(eset)
        print 'old F:',fused
        if 'w' in fused:
            nu = [e for e in eset if not ('w','w') in zip(e,fused)]
            print nu  
            return consistent(nu)
        else: return False
    else: return True

## print 'hello mate'
#badtab = [('[lL][lL]l[[lL]]', (7,8,0,3,0,5,0,1,1,0)), 
#('[[L]][lL][lL][lL]',(9,12,1,4,0,0,6,1,0,0)),
#('[[lL]][L][lL][lL]', (10,11,1,4,0,0,5,1,0,0)),
#('lllll[[lL]]', (5,0,0,1,0,5,0,1,5,0)), 
#('[[Ll]]l[Ll][Ll]', (8,7,0,0,3,0,5,1,1,0)), 
#('[lL][lL][lL][[L]]', (12,9,1,4,0,6,0,1,0,0)),
#('[lL][lL][L][[lL]]', (11,10,1,4,0,5,0,1,0,0)), 
#('[Ll]lll[[Ll]]', (5,5,0,0,2,5,0,1,3,0)), 
#('[Ll][Ll][L][[Ll]]', (11,10,1,1,3,5,0,1,0,0)),
#('[[lL]]lllll', (0,5,0,1,0,0,5,0,5,0)), 
#('l[Ll][Ll][[Ll]]',(9,6,0,0,3,5,0,1,1,0)), 
#('[Ll][Ll][[Ll]]l', (6,9,0,0,3,4,1,0,1,0)),
#('[[lL]]l[lL][lL]', (8,7,0,3,0,0,5,1,1,0)), 
#('[lL]lll[[lL]]', (5,5,0,2,0,5,0,1,3,0)), 
#('[[lL]][lL][lL][L]', (12,9,1,4,0,0,5,1,0,0)),
#('l[[Ll]][Ll][Ll]', (9,6,0,0,3,1,4,1,1,0)), 
#('[Ll][Ll]l[[Ll]]', (7,8,0,0,3,5,0,1,1,0)), 
#('[[L]][Ll][Ll][Ll]', (9,12,1,1,3,0,6,1,0,0)),
#('[Ll][Ll][Ll][[L]]', (12,9,1,1,3,6,0,1,0,0)), 
#('[L][Ll][Ll][[Ll]]', (9,12,1,1,3,5,0,1,0,0)), 
#('[Ll]ll[[Ll]]l', (4,6,0,0,2,4,1,0,3,0)),
#('lllll[[Ll]]', (5,0,0,0,1,5,0,1,5,0)), 
#('[[Ll]]lll[Ll]', (5,5,0,0,2,0,5,1,3,0)),
#('llll[[lL]]l', (4,1,0,1,0,4,1,0,5,0)), 
#('[lL][lL][[lL]]l', (6,9,0,3,0,4,1,0,1,0)), 
#('l[lL][lL][[lL]]', (9,6,0,3,0,5,0,1,1,0)),
#('[lL]ll[[lL]]l', (4,6,0,2,0,4,1,0,3,0)), 
#('llll[[Ll]]l', (4,1,0,0,1,4,1,0,5,0)),
#('[[lL]][lL][lL]l', (6,9,0,3,0,0,5,0,1,0)), 
#('[[Ll]][Ll][Ll][L]', (12,9,1,1,3,0,5,1,0,0)), 
#('[[Ll]][Ll][Ll]l', (6,9,0,0,3,0,5,0,1,0)),
#('[[Ll]]lllll', (0,5,0,0,1,0,5,0,5,0)), 
#('[L][lL][lL][[lL]]', (9,12,1,4,0,5,0,1,0,0)), 
#('l[[lL]][lL][lL]', (9,6,0,3,0,1,4,1,1,0)),
#('[[Ll]][L][Ll][Ll]', (10,11,1,1,3,0,5,1,0,0)), 
#('[[lL]]lll[lL]', (5,5,0,2,0,0,5,1,3,0)) ]

def altERCs(inp):
    if inp != '#lllllll#':
        raise ValueError
    alterc = {} 
    for C,V1 in badtab:
        tempE = []
        for (C2,V2) in badtab:
            tempE.append(ercBit(V1,V2))
        alterc[C] = frozenset(tempE)
    return alterc


#fixtab()